SUBROUTINE make_lattice
USE param
USE global
IMPLICIT NONE
INTEGER :: i,j,k,l,m,n,c1,c2,c3
INTEGER :: n_atom2,foo,cell_no
real(dbl) :: b1,b2,b3

1000 FORMAT(3X,'No. of lattice sites: ',I10,' / No. of atoms: ',I10)

! Lattice info
n_lattice=n_cell*2
ALLOCATE(lattice(n_lattice))
ALLOCATE(atom(n_atom))

OPEN(UNIT=1,FILE=lgfn,STATUS='OLD',ACCESS='APPEND')
WRITE(1,1000) n_lattice,n_atom
PRINT 1000,n_lattice,n_atom
CLOSE(1)

! Set lattice atoms
n_atom2=0
n_cell_plane=nside_int*nside_int
DO k=0,nside_int-1
   DO j=0,nside_int-1
      DO i=0,nside_int-1
         cell_no=n_cell_plane*k+nside_int*j+i
         cell_box_map(i,j,k)=cell_no
         cell_box(cell_no)%n_lattice=0
         cell_box(cell_no)%n_neighbor=0
         DO n=k-1,k+1
            DO m=j-1,j+1
               DO l=i-1,i+1
                  c1=l
                  c2=m
                  c3=n
                  foo=0
                  IF(l.GE.0.AND.l.LT.nside_int) THEN
                     foo=foo+1
                     b1=0
                  ELSE IF(l.LT.0) THEN
                     c1=c1+nside_int
                     b1=-1
                  ELSE
                     c1=c1-nside_int
                     b1=+1
                  END IF
                  IF(m.GE.0.AND.m.LT.nside_int) THEN
                     foo=foo+1
                     b2=0
                  ELSE IF(m.LT.0) THEN
                     c2=c2+nside_int
                     b2=-1
                  ELSE
                     c2=c2-nside_int
                     b2=+1
                  END IF
                  IF(n.GE.0.AND.n.LT.nside_int) THEN
                     foo=foo+1
                     b3=0
                  ELSE IF(n.LT.0) THEN
                     c3=c3+nside_int
                     b3=-1
                  ELSE
                     c3=c3-nside_int
                     b3=+1
                  END IF
                  cell_box(cell_no)%n_neighbor=cell_box(cell_no)%n_neighbor+1
                  cell_box(cell_no)%bx(cell_box(cell_no)%n_neighbor)=b1
                  cell_box(cell_no)%by(cell_box(cell_no)%n_neighbor)=b2
                  cell_box(cell_no)%bz(cell_box(cell_no)%n_neighbor)=b3
                  IF(foo.EQ.3) THEN
!                     cell_box(cell_no)%n_neighbor=cell_box(cell_no)%n_neighbor+1
                     cell_box(cell_no)%neighbor_cell_list(cell_box(cell_no)%n_neighbor)=&
                          n_cell_plane*n+nside_int*m+l
                  ELSE
                     cell_box_map(l,m,n)=n_cell_plane*c3+nside_int*c2+c1
                     cell_box(cell_no)%neighbor_cell_list(cell_box(cell_no)%n_neighbor)=&
                          cell_box_map(l,m,n)
                  END IF
               END DO
            END DO
         END DO
         DO l=1,nbasis
            n_atom2=n_atom2+1
            IF(n_atom2.GT.n_atom) THEN
               STOP 'ERR:make_lattice-n_atom2 cannot be larger than n_atom.'
            END IF
            cell_box(cell_no)%n_lattice=cell_box(cell_no)%n_lattice+1
            cell_box(cell_no)%lattice_list(cell_box(cell_no)%n_lattice)=n_atom2
            ! Initialize lattice info
            lattice(n_atom2)%removed=.FALSE.
!            lattice(n_atom2)%n_near_atom=0
            NULLIFY(lattice(i)%near_atom_list)
            ! Setting lattice position 
            lattice(n_atom2)%px=(i+basis(l,1))*dist_a0_ang
            lattice(n_atom2)%py=(j+basis(l,2))*dist_a0_ang
            lattice(n_atom2)%pz=(k+basis(l,3))*dist_a0_ang
            ! Setting boundary lattice
            DO m=1,n_rlayer
               IF(i.EQ.remove_layer1(m).OR. &
                  j.EQ.remove_layer1(m).OR. &
                  k.EQ.remove_layer1(m).OR. &
                  i.EQ.remove_layer2(m).OR. &
                  j.EQ.remove_layer2(m).OR. &
                  k.EQ.remove_layer2(m)) THEN
                  lattice(n_atom2)%removed=.TRUE. 
                  EXIT
               END IF
            END DO
         END DO
      END DO
   END DO
END DO
! Check
IF(n_atom2.NE.n_atom) THEN
   STOP 'ERR:make_lattice-The number of atoms is not consistent with calculated one.'
END IF
END SUBROUTINE make_lattice
